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Abstract:  A  high-order  locally-corrected  Nystrom  (LCN)  solution  of  a  hybrid  volume/surface  integral  equation 
is  presented  for  the  electromagnetic  scattering  by  complex  targets  that  consist  of  composite  homogeneous  and 
inhomogeneous  materials  and  conducting  objects.  It  is  found  that  for  general  scattering  objects,  the  use  of 
mixed-order  basis  functions  accelerates  the  convergence  of  the  LCN  solution,  eliminates  spurious  charges,  and 
significantly  reduce  the  condition  number  of  the  resulting  linear  system  of  equations.  Fast  iterative  solvers 
based  on  a  multi-level  fast  multipole  method  (MLFMM)  and  the  quadrature  sampled  pre-corrected  FFT 
(QSPCFFT)  are  also  employed  to  enable  the  practical  solution  of  large  complex  targets. 

Keywords:  Electromagnetic  scattering,  integral  equation  methods,  locally  corrected  Nystrom  method, 
fast  solution  methods. 


1.  Formulation 


Consider  the  electromagnetic  interaction  with  a  material  scatterer  composed  of  penetrable  and  conducting 
materials  as  illustrated  in  Fig.  1.  A  surface  separating  volumes  Vj  and  E/  is  denoted  as  St  . .  Let  SG  denote 

the  surface  just  inside  Vh  and  ST  denote  the  surface  just  inside  Vr  Equivalent  current  densities  are  then  placed 
on  all  surfaces  separating  each  material  volume  and  are  defined  as: 

(1) 


J i  i  =hiX-H  + 
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where  hi  and  /L  are  the  unit  normal  directed  into  Vj  and  (A  ,  respectively.  Sine ehj--hi, 

4  Kj=-K=K- 

On  the  surface  of  a  conductor,  only  the  electric  current  density  is  supported.  Consequently, 


J,  „  =  x  H 

l  iP  l  r-+ 


(2) 


(3) 


Each  material  volume  is  assigned  a  background  material  profile (eib,Mib)  ■  I  f( , //, )  ^  (eih,pjb) ,  this  results  in 
the  need  of  an  equivalent  volume  current  density,  expressed  as: 


^Vi  (r)  —  J^ib 


V  Sib 


E(r),  MVi(r)  =  jG)pib 


H(7-) 


Ml)  , 

V  Mib 

where  (sir,pir)  are  the  permittivity  and  permeability  of  the  physical  material  medium. 

If (iy//,)  =  (slb,pih) ,  the  volume  current  densities  are  zero.  Also,  if (sjb,pjb)  =  {sib,pib),  the 


(4) 
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Fig.  1  Equivalent  currents  representing  the  scattered  fields  due  to  a  plane  wave  impinging  on  an  arbitrary 
inhomogeneous  material  scatterer. 


surface  current  densities  on  S,  ,  will  be  zero. 


Each  equivalent  current  density  radiates  in  an  equivalent  homogeneous  material  volume.  For  example, 
J]j ,  M+j,JVi,  Mv,  are  effectively  radiating  in  a  homogeneous  space  with  material  profile  [sib,juib)  ■  Thus,  the 

scattered  field  radiated  by  the  equivalent  currents  in  Vj  can  be  expressed  as: 

Er(jeqMeq)  =  Vib 


Hr  (lq . ^  )  =  K ,(•/„)  +  Tj-%  (Meq ) 


where, 


^(K)  =  -jk, J 
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i  {^eq)  ~  {  V  x  Gt  (r,r')Xeq  (f')dCl' 


(5) 
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where  Q  is  either  a  surface  or  volume,  1  is  the  idem  factor,  or  unit  diad,  G,  ( r,r ')  =  e  JkibV  I  Air  \  r -r'\, 

K  =  ®4etbMtb  and  V*  -  \[mH>  /  £ib  • 

The  proposed  method  of  moment  solution  is  based  on  a  hybrid  surface/volume  integral  formulation  with  the 
appropriate  constraints  applied  on  boundary  surfaces  and  within  material  volumes.  At  a  surface  Sj  .  separating 

two  material  regions,  the  surface  integral  formulation  is  derived  via  a  combined  field  formulation  based  on 


Muller’s  formulation,  leading  to 


I  ■  E" 


+1  ■  £,  E‘] 


1,2]: 
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i sr,j 
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(9) 

(10) 


where  E"W,H'"C  are  radiated  by  impressed  sources  inf],  E":cat , H]“"  are  radiated  by  equivalent  currents  in 
volume  Vi  (similarly  for  subscript  j),  fr  is  the  unit  normal  to  Si  .  and  directed  into  Vi ,  and  7  is  a  test  vector 
tangential  to  5,.  . .  Muller’s  formulation  has  the  advantage  over  the  classical  PMCHWT  (Poggio,  Miller,  Chang, 

Harrington,  Wu  and  Tai)  formulation  [3]  in  that  it  behaves  as  a  second-kind  integral  equation  for  moderate  to 
low  contrast  materials.  Furthermore,  the  hyper-singularity  of  the  L  -operator  in  is  reduced  by  one  order  and  we 
have  found  it  leads  a  reduction  in  the  condition  number  of  the  matrix  arising  from  the  FCN  method. 

On  the  surface  of  a  conductor,  a  combined  field  integral  operator  (CF1E)  is  applied  [4]: 
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where  Zs  is  the  Leonotovich  surface  impedance,  and  t  is  a  test  vector  tangential  to  the  conductor  surface  Si  . 
Finally,  in  regions  where  the  volume  currents  are  non-zero,  the  constraints  applied  are: 

Kc  (?)  =  jr,  {£i(r)lsib  - 1))  -  &r  (?)  (12) 

Hr  (?)  =  MVi  (?)/ (  j ih  ~ !))  “  H-ca'  (?)  •  (13) 

where  ref. 


2.  Discretization 


It  is  assumed  that  the  volume  regions  (as  needed)  and  surfaces  are  discretized  using  curvilinear  cells  that 
represent  all  boundaries  and  material  inhomogeneities  to  high-order  accuracy.  In  this  work,  volume  cells  are 
assumed  to  be  curvilinear  hexahedron  and  surface  cells  are  assumed  to  be  curvilinear  quadrilaterals.  The  current 
density  is  then  expanded  over  each  cell  using  a  set  of  local  vector  basis  functions.  To  date,  LCN  formulations 
have  employed  current  basis  functions  that  are  complete  to  a  polynomial  order  [5,  6].  Here,  the  vector  basis 


used  are  mixed-order  basis  functions,  described  as  [7]: 

jr{u\u2)  =  axPj(u)Pk(u2)ly[g  (j  =  0..p;k  =  0..p-l),  (14) 

J?(u1,u2)  =  a2PJ(u')Pt(u 2)/^r  (j  =  0..p-l;k  =  0..p),  (15) 

on  surface  patches,  or 

J'r  (u,u1,ui)  =  axPj  (u)Pk  (; u1)pi ;  (u^/yfg  (j  =  0 ,.p;k,l  =  0../?-l)  ,  (16) 

jr{ul,u2y)  =  a2Pj(uyk(u2)PI(u2)/^  ( k  =  0..p;j,l  =  0..J3-1) ,  (17) 

j™° (u\u2,u3^  =  a3Pj(ui^Pk(u2^P,(u3}/ *Jg  [j,k  =  0..p-\;  l-0..p )  (18) 


over  a  volume  cell.  Note  that  the  same  basis  functions  are  used  for  magnetic  current  densities.  For  this  choice 
of  basis  functions,  the  charge  density  (derived  from  the  divergence  of  the  surface  current  density)  is  effectively 
described  by  a  polynomial  expansion,  complete  to  order  p.  It  was  shown  in  [7]  by  £ali§kan  and  Peterson  and 
further  demonstrated  by  Gedney,  Zhu,  and  Lu  in  [8]  that  polynomial  complete  current  basis  can  lead  to  spurious 
solutions  for  problems  with  edge  singularities  and  that  this  can  be  remedied  through  the  use  of  the  mixed-order 
basis.  This  is  somewhat  intuitive  for  the  following  reason.  If  one  observes  the  charge  density  at  the  edge 
singularity,  one  expects  the  charge  density  to  exhibit  the  same  leading-order  singularity  as  the  current  flowing 
transverse  to  the  edge  and  as  the  normal  derivative  of  the  current  flowing  normal  to  the  edge.  This  can  be  well 
approximated  by  the  mixed-order  basis.  For  the  polynomial  complete  basis  this  is  approximated  by  the  lower- 
order  terms,  but  there  are  additional  spurious  charge  terms  from  the  highest-order  basis  that  exist.  For  smooth 
scatterers,  this  does  not  occur  and  consequently  either  basis  function  can  be  used  successfully.  However,  it 
appears  that  since  the  charge  density  is  over-specified  by  the  polynomial-complete  basis,  a  poorly  conditioned 
system  results  for  the  EFIE  or  CFIE  operators. 

When  employing  the  mixed-order  basis  functions,  an  appropriate  testing  scheme  must  be  implemented 
within  the  context  of  the  LCN  method.  This  is  done  by  introducing  an  appropriate  quadrature  rule  over  each 
cell.  For  a  surface  curvilinear  quadrilateral  cell,  the  quadrature  rules  are  represented  by  the  product  of  two  one¬ 
dimensional  Guass-quadrature  rules.  For  a  volume  curvilinear  hexahedron,  the  product  of  three  one¬ 
dimensional  Gauss-quadrature  rules  are  used.  For  a  surface  cell,  two  different  quadrature  rules  are  used:  a 
(p  + 1)  x  p  -point  rule  and  a  p  x  (p  + 1)  -point  rule.  Three-quadrature  rules  are  introduced  for  the  volume  cells. 
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Fig.  2.  Curvilinear  quadrilateral  cell  discretization  of  a 
spherical  surface  (24,  5,h-order  cells). 


Fig.  3  Relative  mean  error  in  the  RCS  computed  via  the 
LCN  method  with  ninth-order  basis  and  area  of  a 
sphere  of  radius  koa  —  6  versus  the  interpolation 
order  of  the  24  curvilinear  cells. 


For  each  rule,  a  test  vector  tt  is  introduced.  The  vector  dot  product  of  the  appropriate  integral  operator  is  taken 
with  the  test  vector  tt  at  the  abscissa  points  of  the  quadrature  rule.  This  collocation  scheme  leads  to  a  linear 
system  of  equation.  It  is  noted  that  the  choice  of  the  test  vectors  is  quite  important.  Best  results  are  obtained  if 
the  test  vectors  diagonalize  the  constant  terms  of  the  integral  operators.  Thus,  for  the  surface  integral  operators 
in  (9)  -  (1 1)  suitable  test  vectors  are: 

(sj  —  ~@2’  (s2  ~  ®1’  ^ Pi  —  ®  ’  O' —  5  (19) 

where  a1  is  the  reciprocal  unitary  vector  [9].  Test  vectors  applied  to  the  volume  integral  equations  in  (12)  and 
(13)  are: 

ty=a\  (1  =  1,3).  (20) 


3.  Numerical  Results 

The  above  formulation  has  been  incorporated  into  a  single  program.  The  software  is  capable  of  analyzing 
general  composite  material  scatterers  composed  of  inhomogeneous  or  piecewise  homogeneous  dielectric, 
magnetic  and  conducting  materials  with  an  arbitrary  number  of  regions.  The  surfaces  and  volumes  are 
discretized  with  curvilinear  quadrilateral  and  hexahedral  elements.  The  use  of  high-order  curvilinear  cells  is 
critical  for  the  efficient  implementation  of  a  high-order  solution  scheme.  The  reason  for  this  is  that  a  high-order 
scheme  is  typically  most  efficient  when  applied  to  electrically  large  cells  with  high-order  basis  functions.  The 
reason  for  this  is  that  the  higher  the  order  of  the  basis,  the  steeper  the  slope  of  the  error  curve.  However,  the 
patches  must  also  accurately  represent  the  surfaces  they  are  modeling.  Else,  the  result  can  suffer  from 
discretization  error.  Consequently,  curvilinear  cells  must  be  employed.  To  demonstrate  this,  consider  the 
electromagnetic  scattering  by  a  PEC  sphere  of  radius  a  defined  by  kQa  =  6  .  The  sphere  was  discretized  with  24 
quadrilateral  curvilinear  cells  as  illustrated  in  Fig.  2.  Initially,  the  basis  function  order  was  set  to  p  =  9. 
10x9  and  9x10 -point  Gauss-Legendre  quadrature  rules  were  used  for  the  Nystrom  discretization  of  each 
patch.  Then,  the  order  of  the  cells  was  increased  from  n=  1  to  n  =  12.  The  bistatic  RCS  was  then  computed  for 
the  sphere  and  the  mean  relative  error  was  calculated  relative  to  a  Mie-series  solution  as: 

1 


Mean  Error  =  —  V  - 

N,tt 


!(<SU) 


(21) 
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Fig.  4  Relative  mean  error  of  RCS  for  a  sphere  of  radius  koa  —  6  versus  the  basis  order  p. 

where  Na  is  the  number  of  angles  (360  uniformly  spaced  angles  were  computed).  A  graph  of  the  mean  relative 
error  versus  the  cell  order  n  computed  for  the  MFIE  is  illustrated  in  Fig.  3.  The  error  in  the  area  of  the  sphere  as 
approximated  via  the  Gauss-Legendre  quadrature  rule  is  also  graphed  in  Fig.  2  as  a  comparison.  Initially,  it  is 
observed  that  error  in  the  RCS  and  the  area  follow  the  same  general  trend.  It  also  appears  that  the  minimum 
error  is  reached  when  n  —  p  + 1 .  This  emphasizes  the  need  for  isoperimetric  curvilinear  cells  to  accurately 
represent  the  surface  to  realize  the  full  potential  of  the  high-order  solution  scheme.  The  error  in  the  RCS  was 
also  computed  as  a  function  of  basis  order.  In  each  case,  the  cell  order  n  is  set  so  that  n  =  p+\.  The  RCS  was 
predicted  via  the  MFIE  («  =  0in  (11)),  the  EFIE  (a-  1)  and  the  CFIE  («  =  0.1),  and  the  mean  error  was 
predicted  via  (21).  A  graph  of  the  mean  relative  error  versus  order  p  is  illustrated  in  Fig.  4  for  p  =  2,  4,  6,  and  8. 

Next,  consider  the  EMCC  double  metallic  ogive  [10].  A  288  cell  discretization  of  the  ogive  is  illustrated  in 
Fig.  5.  The  ogive  is  7.5  in  long,  and  has  a  maximum  radius  of  1  in.  The  monstatic  RCS  of  the  ogive  at  9  GHz 
is  illustrated  in  Fig.  6.  This  compares  extremely  well  to  the  measured  and  calculated  data  in  [10]. 

Finally,  consider  the  coated  double  ogive  in  Fig.  7.  The  ogive  is  coated  with  a  thin  dielectric  coating  with  a 
homogeneous  dielectric  of  relative  permittivity  sr  —  9  -  j0.3 .  The  coating  has  a  length  of  8.25  in,  and  a 
maximum  radius  of  1.1  in.  The  RCS  of  the  coated  ogive  at  1.57  GHz  is  illustrated  in  Fig.  8.  Different  orders 
are  illustrated  to  demonstrate  that  the  results  have  converged. 
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Fig.  5  Curvilinear  cell  discretization  of  the  EMCC 
double  ogive. 


Fig.  6  Monstatic  RCS  of  the  conducting  double  ogive 
at  9  GHz. 


Fig.  7  EMCC  double  ogive  with  a  thin  dielectric 
coating. 
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